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The predictability problem in the inverse energy cascade of two-dimensional turbulence is ad- 
dressed by means of high resolution direct numerical simulations. The analysis is done in terms 
of the finite size Lyapunov exponent (FSLE) which is a measure of the growth rate at fixed error 
level. The numerical results are compared with classical closure predictions and good quantitative 
agreement is found. Finally, it is shown that the inertial range predictability properties are not 
affected by the presence of noise induced by the small, not resolved, scales. 
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O ■ I. INTRODUCTION 

<N : 

^ . The knowledge of predictability properties of two-dimensional turbulence is fundamental for estimating the pre- 

^~i ' dictability of atmospheric flow (Leith 1971, Leith and Kraichnan 1972). The characteristic property of two-dimensional 
\^ , flow is enstrophy conservation which force the energy to flow toward large scales. In this inverse cascade regime, dimen- 
sional arguments predict a "5/3" Kolmogorov energy spectrum which is indeed observed in the large-scale atmospheric 
^r^' spectrum (Nastrom, Gage and Jasperson 1984). 

(-^ I From a theoretical point of view, predictability in fully developed turbulence has been investigated as a prototypical 

\^ ' model with many characteristic scales and time. The first attempts to the study of predictability in turbulence dates 
back to the pioneering work of Lorenz (Lorenz 1969) and to Kraichnan and Leith papers (Leith 1971, Leith and 
Kraichnan 1972). On the basis of closure approximations, they were able to obtain quantitative predictions on the 
evolution of the error in different turbulent situations, both in two and three dimensions. Their fundamental papers 
become the backbone for more recent approaches (Lesieur 1997). Because predictability experiments in fully developed 
turbulence are numerically very expensive, to our knowledge there are still no attempt to compare closure results with 
"k^ ■ direct numerical simulations. 

p^ ' In this paper we address the predictability problem for inverse energy cascade of two-dimensional turbulence by 

,-H . means of high resolution direct numerical simulations. It has been recently shown that inverse cascade is not affected 
f^ ' by intermittency corrections and velocity statistics is quasi-Gaussian (Paret and Tabeling 1997, Boffetta, Celani and 
t^^ ,, Vergassola 2000). This makes the problem simpler than the three-dimensional case and we expect that closure-based 
^^ predictions are essentially correct. 

At variance with direct cascade, the inverse energy cascade cannot be observed in decaying turbulence. In order 

to sustain the cascade, a continuous input of energy by random forcing at small scales is necessary. Because we are 

^ ] interested in the intrinsic predictability of the model, we will study the problem assuming known the realization of the 

1—1 forcing. In realistic applications, the forcing represents the small scale dynamics not resolved by the two-dimensional 

^ model (i.e. convectivc motion in the atmosphere). In this case one should take into account also the uncertainty 

^ introduced by the random forcing. Because of the hierarchical structure of the characteristic times, we will see that 

k>( ' ^^^ inertial range predictability properties are not affected by the uncertainty introduced by the forcing. 

rS [ This remainder of the paper is organized as follows: Section || is devoted to a brief summary of the classical closure 
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jrt ' results. In Section 111 and we present the numerical results and their comparison with predictions. Section IV is 



devoted to conclusions. 

II. STATISTICAL TURBULENCE PREDICTABILITY 

Given two realizations of the velocity field u'^^)(x, t) and u'^^)(x, t), a suitable measure for the predictability is the 
error field 
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<5u(x,t)-^(u(2)(x,t)-uW(x,t)) (1) 

from which we define the error energy and the error energy spectrum as (Leith and Kraichnan f 972, Lesieur f 997) 

EA{t) = ^ f\Su{x,t)fd^x^ f EA{k,t)dk. (2) 

Normalization in (|]) ensures that E^ik, t) —> E{k) for uncorrelated fields (i.e. for t -^ cxj). 

The two realizations of the turbulent flow are stationary solutions of two-dimensional Navier Stokes equation in 
which the inverse cascade is sustained by energy injection at small-scale wavenumber kf with a Gaussian random 
forcing. As it is customary in predictability experiments, for most of the results presented here the two realizations 
of the forcing will be the same. In this case the growth of the error will be uniquely due to the deterministic chaotic 
dynamics. At the end of Section II] we will discuss the case in which two different realizations of the forcing are 
implemented. 

Assuming that the initial error can be considered infinitesimal, the magnitude of the difference field starts growing 
exponentially and E^it) ~ -Ea(O) exp(L(2)i) where L{2) is the generalized Lyapunov exponent (in non-intermittent 
systems, L{2) = 2 A, see (Bohr, Jensen, Paladin and Vulpiani 1998). The error growth in this stage is confined at 
the faster scales in the inertial range, corresponding in our situation to the scales close to the forcing wavenumber 
kf, while at larger scales the two fields remain correlated (see Figure nl). At larger times, when EA{kf,t) becomes 
comparable with E{kf), the exponential growth terminates, because the two fields are completely decorrelated at 
small scales. The error growth continues at larger scales in the inertial range, where the two fields are still correlated, 
and the time evolution of the error follows an algebraic law. 

The error transfer towards large scale can be quantified by the characteristic wavenumber kE(t) at which the error 
spectrum is a given fraction of the reference spectrum. Leith and Kraichnan (Leith and Kraichnan 1972), defining 
kE such that the relative error spectrum r{k,t) = E^ik-t)/ E{k) is r{kE) — 0.5, found a self-similar error spectrum 
r{k/kE) in which the time dependence is only through kE{t). 

The dimensional prediction, based on the assumption that the time it takes for the error to induce a complete 
uncertainty at wavenumber k is proportional to the characteristic time at that scale gives, within the Kolmogorov 
scaling 



dt ^ 



(3) 



where B is an adimensional constant. By integration one gets 

fc£(t)-2/3 = k-^/' + Be^'\t ~ to) (4) 

If the inertial range is sufficiently wide, we can assume that ko <C fcs and thus ksit) ~ £~^/^t^'^/^ (Lorenz 1969). At 
wavenumbers smaller than kE the error is still small in comparison with the typical energy, while at larger wavenumbers 
the two fields are completely decorrelated. Thus one can assume that in (0) one has 

Using the Kolmogorov spectrum for E{k) one ends with the prediction 

EA{t)^Get (6) 

The numerical constant G in (o) can be obtained only by repeating the argument more formally within a closure 
framework. In the case of two-dimensional inverse cascade the Test Field model predicts G ~ 4.19 (Leith and 
Kraichnan 1972). 

Under the hypothesis of an initial error localized on infinitesimal scale in a wide inertial range (i.e. A:(0) —^ c»), 
reverting equation (0) is possible to relate the predictability time on a certain scale with the characteristic time taken 
from the energy spectrum: 

T,{k) = is- V3fc-2/3 ^ ^ [^(fc)fc3] -1/2 ^ C^^^f^^ (7) 

where C is Kolmogorov constant (C ~ 6 in two dimension) and T{k) is the spectrum-base characteristic time. This 
result has a practical utility, because in atmospheric turbulence is obviously impossible to perform a "predictability 
experiment" and equation (R) give an estimate for predictability time simply measuring the energy spectrum. 



An alternative approach to predictability problem is the recently introduced Finite Size Lyapunov Exponent analysis 
(Aurell et al. 1996). FSLE is a generalization of the Lyapunov exponent to finite size errors. Within this approach 
one computes the error doubling time Tr{S) which takes for an error of size i5 (with a given norm) to grow of a factor 
r (typically r = 2) From the average doubling time one defines the FSLE as 

It is easy to show that definition (0) reduces to the standard Lyapunov exponent A in the limit of infinitesimal errors 
S ^ (Aurell et al. 1996). As the Lyapunov exponent can be seen like the inverse of the fastest characteristic time 
of a dynamical system, the FSLE can be considered the inverse of the characteristic time at a certain scale. For 
finite errors, the FSLE measures the effective error growth rate at error size S. Let us remark that taking averages 
at fixed time, as in (|6|) is not the same of averaging at fixed error size, as in (|^). This is particularly true in the 
case of intermittent systems, in which strong fluctuations of the error in different realizations can hide scaling laws in 
time. From a numerical point of view, the computation of X{S) is not more expensive than the computation of the 
Lyapunov exponent with a standard algorithm. 

In turbulence predictability, a natural measure of the error is (5 = \/2E/^. The assumption of locality, i.e. that the 
FSLE is proportional to the inverse of the characteristic time at the scale k such that the typical velocity is u{k) ~ S, 
combined with Kolmogorov scaling for the velocities produces the scaling law for the FSLE 

A((5) = AeS-^ (9) 

The constant A relates the energy flux in the cascade to the rate of error growth. Its value is not determined by 
dimensional arguments, but is easy to show that in absence of intermittency and for r ~ 1 it can be related to the 
constant G in (||) by A = ln{r^)/{r'^ — 1)G. In the limit r ^ 1 one gets A ^ G. 

Because of the appearance of the energy dissipation at the first power in (|9|), we expect that this scaling law is 
universal, i.e. not affected by possible intermittency in the velocity statistics (Aurell et al. 1996). The scaling law 
(Q) is valid within the inertial range u{kf) < S < U where u{kf) represents the typical velocity fluctuation at forcing 
wavenumber and U ~ \/2E is the large scale velocity. At large errors J ~ t7, we expect error saturation, i<^A -^ E 
and thus \{5) -^ 0. 

III. NUMERICAL SIMULATIONS AND ANALYSIS 

We have performed extensive direct numerical simulation of the two-dimensional Navier-Stokes equation written 
for the scalar vorticity a;(x, t) —V x u(x, t) — —A'ilj{x, t) as 

dtUJ + J{uj,'ilj)^uALj-aLU + f (10) 

where J represents the Jacobian with the stream function tp. The friction term in (M) removes energy at large scales: 
it is necessary in order to avoid Bose-Einstein condensation on the gravest mode and to obtain a stationary state 
(Smith and Yakhot 1993). Physically this term represents the effect of bottom friction. The random forcing / is 
(5-correlated in time and injects energy at wavenumber kf only. 

Numerical integration of diq) is performed by a standard pseudo-spectral code fully dealiased with second-order 
Runge-Kutta time stepping on a doubly periodic square domain with resolution N = 1024. As it is customary in 
numerical simulations, we use hyperviscous dissipation in order to extend the inertial range. 

Stationary turbulent flow is obtained after a very long simulation starting from a zero initial vorticity field. At 
stationarity one observes a wide inertial range with a well developed Kolmogorov energy spectrum E{k) = Ce^^^k~^^^ 
(see Figure ^) . Structure functions in physical space are found in agreement with the self-similar Kolmogorov theory 
and no intermittency is detected (Boffetta, Celani and Vergassola 2000). 

The perturbed field is obtained from a configuration of the velocity field u^^)(r, 0) in the stationary state, 

u(2)(r,0) = u(i)(r,0) + V2 5u(r,0) (11) 

in which the initial error (5u(r, 0) is very small {E/^{t = 0) ~ 10~^E), and is localized on small scales. The error energy 
spectrum (nh is thus initially localized at wavenumbers greater than the forcing wavenumber. In the inertial range 
the two fields are completely correlated at i = 0. Let us remark that the precise form of the initial error spectrum is 
not important provided that it can be considered infinitesimal, at it will immediately evolve toward the direction of 
the first Lyapunov eigenvector. 



The two configuration are integrated in time according to (|l^) and the evolution of the error 5u{r, t) is computed. 
In this section we present the results obtained using the same realization of random forcing in both the realizations. 
In this way the growth of the error is only induced by the turbulent dynamic. At the end of this section we will 
compare this results with the ones obtained with two independent realizations of the forcing. 

In Figure we plot the time evolution of the error energy (i?A(i)) obtained from direct numerical simulations 
averaged over 20 realizations. The exponential regime is clearly visible at small times, showing that for infinitesimal 
error turbulence behaves exactly as a standard chaotic system. At the very beginning is possible to observe an initial 
recorrelation of the two fields, i.e. a decreasing of E/^(t): it takes a small, but finite time for the initial perturbation 
to align in the direction of the leading Lyapunov exponent. During this time the forcing and dissipation recorrelate 
the two fields. 

Figure shows the computation of the FSLE from our simulations. For very small errors, S < u{kj) corresponding 
to an error spectrum EA{kf,t) << E{kf), we observe the convergence of A(<5) to the leading Lyapunov exponent 
A ~ 1.07. Its value is essentially the inverse of the smallest characteristic time in the system and represents the 
growth rate of the most unstable features. At larger errors 5 > 10"^, we clearly see the transition to the inertial range 
scaling (^. At further larger errors S ~ U ~ 0.1, X{S) falls down to zero in correspondence of error saturation. 

In order to emphasize scaling (O), in Figure we also show the compensation of X{S) with e6~'^. Prediction IM) 
is verified with very high accuracy which allows to determine the value oi A — 3.9 ± 0.1. With the present value 
of r ~ 1.057, this corresponds to a value G = 4.1 ± 0.1. The physical picture we obtain is that the creation of 
uncorrelated energy in the inertial range due to chaotic dynamics is about 4 times faster than the energy transfer 
rate. Our numerical result is in remarkable agreement with the old prediction obtained within the Test Field Model 
closure which gives G = 4.19 (Leith and Kraichnan 1972). 

The physical meaning of G is the ratio of the rate of uncorrelated energy production to the rate of energy injected 
by the forcing and transfered to large scales e. The fact that G > 1 shows that the uncorrelated energy is not simply 
"transported" with the energy through the cascade, but there is an effective production of error at each scale, due 
to the chaotic dynamic. The constant rate of error-energy growth is not an immediate consequence of the constant 
energy fiux in the inverse cascade, but is the effect of the dimensional hypothesis that the times for energy transfer 
and error growing at a fixed scale should follow the same scaling law. Moreover G > 1 suggests that (0) should be 
unchanged when one considers two independent realizations of the forcing, as we will see below. 

The time evolution of the error energy spectrum (Figure m) shows a sharp front of error thatpropagates from small 
to large scales justifying a posteriori the assumption (^. From the error spectra in Figure yj one can compute the 
characteristic wavenumber ksit) defined as r{kE) = 0.5. On scales smaller than ksit)^^ the two fields are already 
decorrelated, while on larger scales the error can still be considered small, so we can define implicitly the predictability 
time Tp(k) at certain wavenumber k from kslTp) — k. 

Figure |j shows the evolution of kE{t)~'^^^ compared with the best fit based on (|j). The result is rather noisy 
because of the discrete character of the wavenumbers, nevertheless it is possible to estimate the value of the constant 
in (|j) _B = 0.43 ± 0.02. In equation (0) B relates the predictability time to the characteristic time based on the 
spectrum. From our simulation we have Tp{k) ~ 5.7t(A;). 

In order to check the prediction of the self-similar error growth we rescale all the error spectra with the energy 
spectrum, obtaining the relative error spectra r{k,t) = EA{k,t)/E{k). Rescaling the wavenumbers with the charac- 
teristic wavenumber kE{t) the error spectra collapse in a similarity error spectrum r{k/kE) (Figure H). This means 
that in order to characterize statistically the error growth one needs to know only the evolution of fc£;(i). Closure 
computations also predict that the shape of the error spectrum in the limit fc ^ is EA{k,t) ~ k^^ that for the 
similarity error spectrum gives r{x) ^ x^"^'^ which is indeed observed in our simulations. 

A. Different forcing 

All the results presented above are obtained forcing identically the two realizations of the turbulent field. In principle 
this is not consistent with a realistic application like modeling large scale atmospheric flow because in that case the 
forcing is due to the motion on small scales that are non resolved. In order to simulate this continuous injection of 
error it is more correct to use two independent realizations of the random forcing. To check in which way this can 
affect the results discussed above, we have repeated the same numerical experiment with two independent forcing. 

The main difference is that in this case the initial exponential growth of the infinitesimal error is non present. 
When the scales smaller than the forcing scale are non completely decorrelated the main source of error is due to 
the different forcing, not to the chaotic dynamic. If there were no dynamical effects we could expect to observe a 
"diffusive" behavior and the error energy should grow linearly as in (|6|) but because the source of the error is the 
energy input one should have G ~ 1. This regime is indeed observed (Figure 0) but is worthwhile to remember that 



it is physically of little interest: we are using different forcing to mimic the unresolved small scale motion and so it is 
natural to assume that the two fields are completely decorrelated at those scales. 

When this condition is satisfied, and the error is localized in the inertial range, it is not possible to distinguish the 
case with a different forcing from the one with one forcing. From the analysis of the FSLE (Figure g) it is evident 
that in this range the scaling law (H) is verified with the same value of the constant A. All the inertial range results 
about the self similar growth of error are recovered. This demonstrate that the continuous injection of error due to 
an unresolved small-scale motion do not produce effects when the two turbulent fields are already decorrelated on 
such scales. This is a consequence of locality: the error growth in the inertial range is completely determined by the 
non linear dynamic, and it is not sensible on what happens on smaller scale. 

IV. CONCLUSION 

We have studied the predictability problem in the inverse cascade of two-dimensional turbulence. Using the Finite 
Size Lyapunov Exponent we have shown that after an initial exponential growth, a linear growth of the decorrelated 
energy sets in. This decorrelated energy flows to large scales in self similar way through an inverse cascade with a 
constant flux that is about four times faster than the energy flux. 

This behavior do not change also adding other possible perturbation on small scale, thus the results obtained in a 
numerical experiment with initial condition uncertainty can be applied to a wider class of problems. As an example 
the predictability time of mesoscale in the atmosphere can be estimated easily with a simple measure of the energy 
spectrum. 

All the results presented are in strong agreement with the prediction of closure. At least from the point of view 
of predictability, two-dimensional turbulence thus seems to be very well captured by low-order closure scheme. As 
a consequence we can exclude, on the basis of our numerical findings, the existence of intermittency effects in the 
inverse cascade of error. This is a result which is probably of more general interest than the specific problem discussed 
in this article. 
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FIG. 1. Stationary energy spectrum E{k) (thick fine) and error spectrum E^ik, t) at time t — 4.6, 5.5, 7.1, 10.0, 15.6. fc/ = 320 
is the forcing wavenumber. In the inset we plot the compensated spectrum £~'^'^k^'^E{k). 



FIG. 2. Error energy {-EA(i)) growth averaged over 18 runs. Dashed line represents closure prediction (pi), dotted line is the 
saturation value E. The initial exponential growth is emphasized by the lin-log plot in the inset where the initial decreasing of 
the error is also observable. 

FIG. 3. Finite size Lyapunov exponent A(5) as a function of velocity uncertainty 5. In the simulation with identical forcing 
for the two fields (+) the asymptotic value for 5 —* Q gives the leading Lyapunov exponent of the turbulent flow. In the case 
with different forcing (x) the infinitesimal regime is unphysical. In the inertial range, 5 > 10"'^, the behavior of \{S) is identical 
for both case. Dashed line represent the prediction (bl). In the inset we show the compensated plot \{5)S^/£. The line represent 
the fit to the constant A ~ 3.9. 

FIG. 4. Time evolution of the characteristic wavenumber fc_E(i)~^" compared with the dimensional prediction 
fcE(i)~^ = kg + Be^'^{t — to) where to ~ 6.0 is the time that the error takes to reach the inertial range at fco — 133.0 

FIG. 5. Similarity error spectrum r{k/kE) = EA{k/kE,t)/ Eik/ks)- All the error spectra rescaled collapse together within 
the error bands. The log-log inset shows the asymptotic behavior r{x) ^ a;^**" for x — » 
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